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Hypotheses that do not have the standard bilinear form 0 = CBM = 0 occur naturally in the analysis of 
repeated measurement designs. An expanded class of tests of the form y/ = Tr(GQ) = Q, called extended linear 
hypotheses, provide a richer class of parametric functions. In this paper we show how to analyze double multivariate 
(DMM) and mixed multivariate models (MMM) using SAS. We extend the analysis to extended hypotheses and 
develop a new approximate test of extended linear hypotheses for MMM designs that does not require multivariate 
sphericity, only a general Kronecker structure. 
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1. Introduction 



In the analysis of repeated measurements, one may be interested in modeling changes in the matrix of means, 
testing hypotheses regarding the location parameters, and modeling the covariances among the repeated 
measurements. Standard statistical analysis of such designs under multivariate normality use the multivariate 
analysis of variance (MANOVA) model or the general growth curve model, the GMANOVA model developed by 
Potthoff and Roy (1964), Timm (1980) and von Rosen (1991). These models may not be used to test all hypotheses 
about the parameter matrix B of means since certain linear hypotheses about B may not have the simple bilinear form 

//:0 = CBM = O. (1.1) 

For example, if 

<^ 2 ! M22; 

for a simple 2x2 repeated measurement design, the hypothesis 

=0and/t2! “Mi2 =0 (1-2) 

may not be put in the form of (1.1). To solve this problem, Mudholkar, Davidson and Subbaiah (1974) introduced 
the extended linear hypothesis 

//r:rr(G0) = O (1.3) 

for all GeT where F is a well defined set of pxg matrices and Tr{ ) denotes the trace of a matrix. More recently, 
Hecker (1987) introduced the completely general MANOVA model (CGMANOVA) to analyze hypotheses like 
(1.2). In this article, we use the full rank model to analyze hypotheses involving means for the double multivariate 
model (DMM) and the multivariate mixed model (MMM) using SAS (1990) software, analyze hypotheses of the 
form (1.3) for a repeated measurement design, and develop a new approximate test for extended multivariate MMM 
hypotheses without assuming multivariate sphericity. 

2. Repeated Measurements Designs 

. Models for the analysis of repeated measurements are wide and varied. For a comprehensive overview, see the 
recent works of Diggle, Liang and Zeger (1994), Lindsey (1993), and Longford (1993). In this section we review the 
standard repeated measures design under multivariate normality where a vector of p-variates is observed over q 
occasions discussed by Boik (1988, 1991). 

For a multivariate repeated measures design, we observe p-variates over q conditions for N subjects where 
/= 1,2, = 1,2, ...,/ 7 , and/: = 1,2, ...,^. The data matrix Y of order for the design is displayed in 
Figure 1,1 
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Figure 1.1 Data layout for multivariate repeated measures design. 

Using matrix notation, the double multivariate linear model (DMM) for the design may be expressed: 

Y = XB + U (2.1) 

where is the data matrix, is the full rank design matrix, rank (X) = k, is the matrix of fixed 

unknown location parameters, and U^^^ is a matrix of random errors. Each row of U is assumed to be 
independently, normally distributed with common covariance matrix ^pqxpq • Using the standard vec( ) operator that 
stacks columns of a matrix, we write that 

vec(U)--A^[O,(Z0I^)] (2.2) 

where 0 denotes the standard, right Kronecker matrix product (A 0 B = {a,yB}), Graham (1981). 

Given (2.1) and (2.2), we are primarily interested in estimating B and testing hypotheses of the general form: 

//:0 = CB(Ip0A) = O (2.3) 

where is the hypothesis test matrix of full row rank, rank (C) = g, and ^gxu ^s a within subjects matrix of full 
column rank u, rank (A) = u. Without loss of generality, A is assumed to be semi-orthogonal, A 'A = 1^. 

The least squares (LS) and maximum likelihood (ML) estimator of © is 

© = CB(Ip 0 A) where B = (X’ X)‘‘ X'Y. (2.4) 

The estimator © is unbiased and the covariances of the elements are given by 

cov( vec 0) = n 0 C(X'X)“' C' (2.5) 

where SI = (1^ 0 A')2(Ip 0 A). Given (2.2), the 

vec (0) ~ A^(vec 0, ft 0 C(X'X)“' C') (2.6) 

SO that the distribution of © depend on Z through SI. When SI has special structure, mixed model procedures may 
be used to test (2.3). 

When SI is known, the likelihood ratio (LR) statistic for testing (2.3) is 

= vec'(©)[fl 0 C(X'X)“' C']“' vec (0) (2.7) 

and vec'(0) = (vec ©)'. X^ has a chi-square distribution with noncentrality parameter S : 

~ x\pgu,S) 

S = vec'(0)[fl 0 C(X'X)“' C']“' vec (0). 



( 2 . 8 ) 



(2.9) 



However, is usually unknown. Then standard tests of (2.3) depend on the characteristic roots of HE”' where 

E = (Ip ® A)' Y'(Ijv - X(X'X)-' X')Y(Ip ® A) 

H = 0'(C(X'X)”'C')”'0 

are the error and hypothesis sum of squares and cross-products matrices. An unbiased estimator of 12 is 
Q = E/(_N — k). Given (2.2), E and H have independent Wishart distributions: 

E-IV (A^-I:,aO) 

, ( 2 . 10 ) 

H - lVp„(g,fi,f2-''P) 

with noncentrality matrix 4^ = 0'(C(X'X)”'C')”'©. The standard multivariate test statistics: Wilks’s A, Roy’s 
largest root, the Bartlett-Lawley-Hotelling trace and the Bartlett-Nanda-Pillai trace criteria all depend on the 
characteristic roots of HE"^ The Bartlett-Lawley-Hotelling trace statistic, also called Hotelling’s generalized 
statistic, has a simple form. Using the fact that the vec(ABC) = (C'0A)vec B and the Tr(AB)= vec'(A')vecB 
implies that the rr(AZ'BZC) = 7r(Z'BZCA) = vec'(Z)(A'C'0 B)vec(Z), it is easy to show that 

r? = (A^- )k)rr(E”'HE”'H) = (A^- )k)ZA, 

<• ( 2 . 11 ) 

= vec'(0)(Q®C(X'X)”' C')”' vec (0) 
where A,- are the characteristic roots of |H — AE| = 0 or of HE'' . 

Pillai (1960) tabled percentage points of U ='L^i = T„ ! (N—k). Alternatively, one may also employ an F- 

i 

distribution approximation, SAS Institute (1990, p. 19) developed by Pillai and Young (1971). The statistic is 

2(£rt + l^^ - F[j(2m+j+l),2(jrt+l)] (2.12) 

s^(2m+ s+ 1) 

where s= min(g,/?M), m = (lg-/?Ml-l)/2, and n = {N-k- pu-l)/ 2 has an approximate F-distribution when H in 
(2.3) is true. For large departures from H, appears to be best, this is not the case for small departures, Seber 
(1984, p. 415). Power for each of the tests for specific alternatives is easily approximated using a SAS IML routine 
developed by Muller, LaVange, Ramey and Ramey (1992). 

If the matrix SI satisfies the multivariate sphericity assumption: 

n = = (2.13) 

where is the covariance matrix of the vector of p-variates, Boik (1988) showed that (2.13) is the necessary and 
sufficient condition for obtaining exact MMM tests from the DMM. Thomas (1983) developed a test for 
multivariate sphericity. Boik (1988) improved the accuracy of the test by using Box’s correction factor. 

Boik (1991) using the generalized trace operation r^(W) s {rr(W^y)} where Wpuxpu = is a wXw 

submatrix of W developed by Thompson (1973), shows how one may obtain MMM tests from the DMM. The tests 
depend on the characteristic roots of 

H* = r^(H) and E* = 7^(E). (2.14) 

For additional detail, see Boik (1988) and Timm (1980). 




3 



(2.15) 



One can also develop a general MMM directly from (2.1). The model for the subject is 

y, =(Ip^<8>x;)vec(B) + Ui 

where x- and u,' are row vectors of X and U, respectively. To construct a mixed model, we let 

u,. =Z,A,. +e,. 



where Z,- is ap^x /i matrix of covariates of rank /i. A, and e, are independently normally distributed: 

A, ~ ^,(0,D,) 
e, ~;V^^[O,(Z,0I^)] 

where and are positive definite matrices. Then from (2.2) and (2.16), we have that 

Yi ~ A^(Ip^<8>x')vec(B),ZJ 
Z,=Z,D,Z' + (Z,0/,). 

Reinsel (1982) in the development of a multivariate growth curve model with a multivariate random-effects 
covariance structure for subjects proposed a model similar to (2.17). To obtain Reinsel’s model, we associate 
Zf s 0 so that h = p and 



(2.16) 



(2.17) 



Z,=Z = (Dp0O;) + (Z,0IJ 



(2.18) 



where is a pxp positive definite matrix. The structure of Reinsel’s matrix is more commonly called multivariate 

compound symmetry, a sufficient but not necessary condition for exact MMM tests. When Z satisfies (2.13), we 
have Scheff6’s MMM which also requires homogeneity of the covariance matrices. Given (2.13), one may directly 
test (2.3) using (2.17), the mixed MANOVA model with homogeneity of the Z,. 

Rearranging rows and columns in (2.18), observe that (2.18) is a special case of the more general linear form 

Z=Gj0Zi + G2 0Z2 (2.19) 

where Gj and G 2 commute with each other. Tests for matrices that have the linear structure given in (2.19) where 
the matrices G, are known and may or may not commute are reviewed by Krishnaiah and Lee (1980). With the rows 
and columns rearranged. Cl has the form 

Q = (A0I^)'Z(A0I^) (2.20) 



so that the sphericity condition becomes 

Q = I„0Z,. (2.21) 

Given the general form (2.19) with Cl defined in (2.20), we see that Cl may be factored: 

Q = Z„0Z, (2.22) 

where Z„ and Z^ are positive definite matrices and Q does not satisfy (2.21), the multivariate sphericity condition. 
Given that Cl can be factored and estimated, we develop a large sample test of (2.3) given (2.22) as an alternative to 
^-adjusted MMM tests which on the whole are less efficient than DMM tests, Boik (1991). Ameniya (1994) has 
developed a test procedure for testing that of subset of Cl is zero, using the Bartlett-Lawley-Hotelling criterion. 



3. Extended Linear Hypotheses 

Not all multivariate hypotheses have the general bilinear form given in (2.1). More generally, one wants to 
consider all parametric functions v^ = d'0f for arbitrary vectors d'andf. However, = rr(f d'©)= rr(G 0) 
where the rank (G) = l. This motivated Mudholkar et al. (1974) to consider the more general decomposition of 
( 2 . 1 ): 

Wi-: n {W(G):rr(G0) = O} (3.1) 

Ger 

where F is the set of all mxg matrices defined: 



r = {G:G = Z A, G,,A/ s real}. 

1=1 

The hypotheses in (3.1) are called the extended linear hypotheses of order v for known matrices = 1,2, ... , v. 

To test (3.1) when v = 1, Mudholkar et al. (1974) developed a union-intersection (UI) test of 

//:V^ = 7r(G0) = O (3.2) 

using an argument involving symmetric gauge functions, Krishnaiah, Mudholkar and Subbaiah (1980). Their test 
statistic is 

t/(G) = rr(G0)/{rr(GC(X'X)”'C'G'E)}''2 (3.3) 

where E is defined in (2.9), and C is given in (2.3) for the overall test. Letting be the critical value for the F- 
approximation in (2.12), 1-a simultaneous confidence intervals for y/ are: 

y/eiif± c4rr(GC(X'X)-'C'G'E)}''^ (3.4) 

^ * 

for all matrices G„^g ¥ = Fr(G0). 

To test any intermediate extended linear hypothesis requires finding the supremum in (3.3). Given 

V 

G|,G 2 , ... ,Gy where G = Z G, for some vector A = (A,) it is easily shown that 

1=1 

{t/(G)}^ = sup{t/(G)}^ = T T"' T (3.5) 

A 



(3.6) 



where the elements of t and T are: 

T, = rr(G,©) 

= rr(G,.C(X'X)-'C'G;.E) 

Mudholkar et al. (1974). Hence, the intermediate hypothesis is rejected if {t/(G)}^ exceeds the square of the critical 
value for the test of H in (3.2). 

To illustrate several extended linear hypotheses, consider a one-way MANOVA design involving three groups 
and three repeated measures so that the parameter matrix is 



Mi2 



B = 



M21 M22 M23 
<^31 M32 

where the primary hypothesis of interest is the equality of group means: 



(3.7) 



5 



7 



H, 



G*' 



^ 11 ^ 
Mil 




^21^ 




f II ^ 

M 31 


Mi2 


= 


M 22 


= 


M32 


<Mi3> 




<M23; 




^M33; 



To test H . , one may select C = and A = A^to test the overall hypothesis where for example, 



C = 



O ‘ 

^1 -1 



0 



and A„ = I3 



(3.8) 



.0 1 - 1 ; 

SO that we may represent//^, in the general form C^B = 0. Upon rejecting , suppose we are interested in 
the contrast 

V “‘■^22 “‘■^33 “(Mi2 “•"M2i)^^“(Mi3 “•"M3i) ^ ^ ”(^23 “•‘^ 32 ) (^*9) 
which compares the diagonal means with the average of the off diagonals, a situation that may be of interest in a 
repeated measures design. Contrast (3.9) may not be expressed as a bilinear form, y/ = c'Ba, Letting 

^ 1 . 5 ^ 

G= -.5 .5 
,-.5 1, 

observe that y in (3.9) may be expressed as an extended linear hypothesis since 

Y = 7r(G 0) = 7r(GC,BA J = 7r(A,GC,B) = 7r(G*B) 
where and A^ are defined in (3.8) and G* = A^GC^. Thus, in terms of the original parameters, we have that 

f 1 -.5 -.5^ 



G* = 



-.5 1 -.5 



[-.5 -.5 



1 



where each of the coefficients in each row and column of G* sum to one. A contrast of this type involving rows and 
columns in a standard two-way design is called a generalized contrast by Bradu and Gabriel (1974) and Boik (1993). 

Testing that the contrast in (3.9) is zero is achieved using a single matrix G. Alternatively, suppose we are 
interested in the following hypothesis 

Mil = M21 

H: /i,2 = H 22 = M32 (310) 

M23 = M33 

which is not a general linear hypothesis of the form CBM = 0, but an extended linear hypothesis. To see this, 
observe that hypothesis (3.10) is equivalent to testing 

Tr{G 0) = rr(G,. C, B A, ) = rr(A,G,C,B) = Tr(G,^ B) = 0 

for i = 1, 2, 3, 4, where 





fl 0 ^ 




^0 0 ^ 




^0 0^ 




"0 0^ 


Gj = 




0 0 

^0 0 


,G 2 = 


1 0 

,0 0 , 


»^3 = 


0 1 
,0 0, 


,G 4 = 


0 0 

<0 1. 



To test (3.10) requires obtaining the supremum in (3.3) as outlined in (3.5). 
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The examples we have illustrated have assumed A„ = I. We now consider the test of parallelism for a repeated 
measures design so that we have 



1 -1 0 ^ 
0 1 -1 



Vii 1 0 



-1 1 
0 -1 



= 0 



^^2\ A^23 

VA^31 A^32 P-'i'ij 

C„BA„= 0 . 

Following the overall test, suppose we were interested in the sum of the following tetrads; 

¥ = (f^2\ +Mn “ A^22) + (A^32 + P-22~ i^2l)- 

Such a contrast may not be expressed using the simple bilinear form, \(f = c'Ba . However, selecting 

fO O 



( 3 . 11 ) 



G = 






we see that the contrast 



where 



Xjf = MG©) = MGC^BAJ = MA^GC^B) = MG^B) 






G* = 



0 1 -1 

1 -2 1 

-1 1 0 



is a generalized contrast matrix. 

As our last illustration, we let C„ = I and for a repeated measures design suppose we are testing for equality of 
conditions: 



H,.: 



Expressing the hypothesis in the linear form, we have that C„ = I3 and 

f 1 0 ^ 

-1 1 



fu y 

Mil 




I'll ^ 

y-\2 




f II \ 

y-\3 


f^21 


= 


f^22 


= 


f^23 


<^31 > 








^^33 > 



K = 



0 -1 



Following the overall test, one may be interested in the contrast 

¥ = {l^n-l^n) + {1^22-1^23) + (^33)- 

This contrast is again an extended linear hypothesis and tested with 



G = 



( 3 . 12 ) 







V\ 




( 1 


0 


1> 


(1 


0 










or G* = 


-1 


1 


0 


U 


1 


1 










V 




/ 


. 0 


-1 


-1> 



in the expression y/ = Tr(G©) = Tr(GC^BA^) = Tr{G * B). 
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4. SAS Illustrations 

DMM/MMM - Example 

To illustrate the application of the DMM and the MMM, we use the dental data analyzed by Timm (1980), 
Thomas (1983), and Boik (1988, 1991) using SAS. The data were obtained from Dr. Thomas Zullo in the School of 
Dental Medicine at the University of Pittsburgh. The data matrix Y consists of nine subjects randomly assigned to 
two orthopedic treatments. At each of three activator treatment conditions (occasions), three dependent variables 
were observed representing three measurements associated with the adjustment of the mandible. 

The parameter matrix B has the matrix form 

1 2 
12 3 12 3 



1 



3 

2 



variable 

conditions 



(4.1) 



B = 



l^n l^\3 



/^17 /^18 f^\9 

1^27 /^28 1^29 J 



/^14 /ti5 

V /^21 /^22 /^23 1^24 /^25 /^26 

The first hypothesis of interest for the DMM is whether the profiles for the two groups are parallel. Is there an 
interaction between treatment groups and conditions? The hypothesis is 

HgC- /il7-/il9.Ail8-Ail9) 

= {^^2l ~P^23’ A^22“/^23» ■■■ » “A^29» ~A^29)- 



(4.2) 



The matrices C and A are 



C = (l,-1) A = 



1 



0 

1-1 -1 



(4.3) 



in the general form CB(I® A) = CBM = 0. Normalizing A so that A'A = I and nesting differences in conditions 
for each variable, we may write M: 

0 0 0 0 0 0 

0 0 0 0 0 0 

.707 0 -.707 0 0 0 

-.408 .816 -.408 0 0 0 

0 0 0 .707 0 -.707 

0 0 0 



M' = 



.707 0 -.707 

-.408 .816 -.408 

0 0 0 

0 0 0 

0 0 0 

0 0 0 



-.408 .816 -.408 



(4.4) 



so that the normalized first difference contrasts are followed by the normalized second difference contrasts. To test 
for vector differences in conditions the hypothesis becomes 



Hr 



f // > 
Mil 










^^ 2 \ 




^^22 




^23 






^^\5 






^^ 2 A 




^^25 




f ^26 






^12 




f ^\9 


^^27 > 




^A^28> 




^^29 > 



(4.5) 
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To test (4.5), the hypothesis matrices are: 



C = I 2 and A = 



1 0^ 
-1 1 
0 -1 



(4.6) 



To test for differences in Groups, // . , the overall hypothesis is 



Hq,: fii = fi2 (4.7) 

and the hypothesis test matrices are 

C = ( 1 ,- 1 ) and A = I 3 (4.8) 

Tests (4.5) and (4.6) do not require parallelism of profiles. Given parallelism, tests for differences in conditions 
and groups are written by averaging over the other factor: 



Hr.: 



( 3 
7=4 






2 /^ 2 y /3 
7=4 



Z fly/ 3 

KM 



(4.9) 



He- 



2 

r/2,1/2 




/ 2 > 

2/2,1/ 2 




/ 2 

2/2,3/ 2 

f— ,1 


1—1 

.i/2,.4/2 

1—1 


- 


1=1 

2/2,5/ 2 
1—1 


- 


1=1 

2/2,6/ 2 
1^1 


2/2,7/ 2 
Vi=l ) 




2/2,8/ 2 

U‘=l ) 




2/2,9/ 2 
^i=i y 



which simplify the interpretation of the tests. The corresponding hypothesis test matrices are: 

"1/3^ 

Cc = (1,-1) and A = 



1/3 

Vl/3y 



Cc = (l/2,-l/2) and A = 



1 0 

0 1 

-1 -1 



(4.10) 



Normalization of the post matrix M associated with (4.10) permits one to obtain the MMM from the multivariate 
tests given parallelism and multivariate sphericity, Timm (1980) and Boik (1988). This is not the case for the tests 
and//^. 

Utilizing SAS, the code for these hypotheses is provided in Appendix A, Program l.sas. Reviewing the output 
for the SAS run, one may construct Table 1 using the Bartlett-Lawley-Hotelling trace criterion and associated F- 
ratios. 
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Table 1 
DMM OutDut. 


Zullo Data 




Hypotheses 


[/ = 7’//v. 


df 


GC 


0.7153 


(1.2, 4.5) 


G* 


1.3684 


(1,3.5, 3) 


C* 


29.3173 


(2, 1.5, 4.5) 


G 


0.1314 


(1,0.5, 6) 


C 


28.6021 


(1,2, 4.5) 



F 


df 


p-value 


1.3114 


(6, 11) 


0.3292 


1.2163 


(9, 8) 


0.3965 


24.4311 


(12, 20) 


<0.0001 


0.6132 


(3, 14) 


0.6176 


52.4371 


(6, 11) 


<0.0001 



Because the test of parallelism is not significant, tests of conditions and groups given parallelism are appropriate. 
From these tests, one may obtain the MMM tests of Hq and Hqq which are valid given multivariate sphericity. The 
MMM test of group differences given parallelism is identical to the DMM since the test of Hq does not depend on 
the sphericity assumption. 

To test for multivariate sphericity, we may use the procedure developed by Thomas (1983). Thomas showed 
that the LR statistic is 

A=IE1^'^/Im“’ (4.11) 
1=1 

where w = rank(A) or w = rank(M)/p. The asymptotic null distribution of —2 In A is a central chi-square 
distribution with degrees of freedom 

/ = p(«-l)(pw + p+l)/2. (4.12) 

u 

Letting a,(i = 1,2, ... ,p«) be the eigenvalues of E, and )3,(i = 1,2, ... ,p)the eigenvalues of Z E,-,-, a simple form of 



(4.11) is 



p pu 

W = -2\nX = N[u'L ln(i3,)- Z ln(a,)]. 
1=1 1 = 1 



(4.13) 



When p and/or u are large relative to N, the asymptotic distribution of W may be poor. To correct for this, Boik 
(1988) using Box’s correction factor for the distribution of W showed that the 



P(W^W„)=Pip*W^p*W„) 

= (1 - co)P(Xj ^p*W„)+ Q)P(Xj^^ ^p*W„) + 0(vf ) 



(4.14) 



where 

p =l-p [2p\u^ - 1) + 3p(u^ - 1) - («2 - 1)] / I2ujv, 

p*=pvJN (4.15) 



\lipu-l)pu(pu+l)(pu + 2) 




(p-l)p(p+l)(p+2) 




■/(I-P)'ll 


[[ 2Av] 




24m^v^ 




^ Jj 



V, = N- rank(X), and / is defined in (4.12). Hence, the p-value for the test of multivariate sphericity using Box’s 
(4.16) correction becomes 
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Pip*W^WJ = {l- o))P{Xj ^W„) + 0)P(xj^4 S ) + 0(v;^ ). 

Using (4.14) and (4.15), the IML code is given in Program l.sas for the test of multivariate sphericity. From 
the program output, 

W = -2 In A = 74.367228 with t//=15 

and the p-value for the test is 7.365 x 10’*®. With Box’s correction, p*W=54.742543 with p- value of 2.7772 x 10'^ . 
Hence, the MMM assumption is not satisfied with Zullo’s data. Thus, one should use the DMM for these data. 

Continuing, assuming the test of multivariate sphericity was not significant, we show how one may analyze 
the DMM data using the MMM. In addition, if one could not test for sphericity, we illustrate how one may construct 
the e-adjusted MMM tests statistics developed by Boik (1988, 1991) which extends the Greenhouse-Geisser type 
procedure to the MMM. While one may directly derive the MMM test matrices from the normalized DMM tests, as 
illustrated in Timm (1980), it is more convenient to use the SAS PROC GLM to obtain the test statistics and 
hypothesis test matrices. 

For the MMM analysis, the data must be reorganized as a mixed model. Subject is a random factor nested 
with the fixed group factor, and groups are crossed with the fixed condition factor. To perform the tests in SAS, we 
must use the RANDOM command to calculate the expected mean squares to locate the appropriate error sum of 
squares matrices for the multivariate tests. From the expected mean square output, we see that the subject within 
group effect is the error term to test for group differences given parallelism. The commands for the analysis become: 
model y l-y3 = group subj(group) cond cond*group; 
random subj (group); 

contrast ‘Group’ group l-l/e=subj (group); 
manova h = cond group*cond/printh printe; 

as illustrated in Program l.sas. The option e on the contrast statement defines the error matrix for the test of groups. 
The “overall” error matrix is being used to test for conditions given parallelism and the test for parallel profiles, the 
condition by group interaction. The MANOVA table for the multivariate mixed model analysis is shown in Table 2. 



Table 2 

MMM Analysis. Zullo Data 



Hypothesis 




df 


F 


df 


p-value 


Group/Parall 


0.1314 


(1,0.5, 6) 


0.6132 


(3, 14) 


0.6176 


Cond/Parall 


13.7514 


(2,0, 14) 


66.4651 


(6, 58) 


<0.0001 


Group X Cond 


0.1907 


(2,0, 14) 


0.9218 


(6, 58) 


0.4862 



That the multivariate mixed model results may be obtained from the DMM, observe that for the test of 
conditions given parallelism that the hypothesis test matrix 
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H* = 



^152.93 iSym)> 

95.79 62.73 
^ 14.75 856 159 ^ 



is the average of Hu and H 22 where 





"148.02 




Symf 




" 4.89 


Sytn' 


H„ = 


96.32 


62.67 




«22 = 


-053 


0.06 




^ 13.38 


8.71 


1.21 j 




, 1.36 


-0.15 0.39 ; 



are the submatrices of H for the test of given parallelism. The degrees of freedom for the mixed model analyses 
are vJJ = rank (M) / /? = rank(A) =1-2 = 2 and v* = rank(M) / p = rank(A) = 16 • 6 / 3 = 32. The 

results for the test of cond x group interaction follow similarly, Timm (1980) and Boik (1988). From the MMM 
results, one can also obtain the univariate split-plot univariate F-ratios, a variable at a time. 

When multivariate sphericity is not satisfied, one cannot perform a MMM analysis. Boik (1988, 1991) 
developed an e-adjustment for the tests of conditions and conditions x group interaction following the work of Box. 
Because the multivariate tests may be approximated by a Wishart distribution when the null hypothesis is true: 

E* ^ Wpiv* = a/v„aO) 

H* ~ Wp(v„* = £«v„n-'A,0) 



Boik showed how one may use an estimate of B to adjust the approximate F tests for the multivariate criteria where 

v2l r ^ xn2 



Tr 



e = 



[I®*] 






(4.17) 



M = rank(A) and andvg are the hypothesis and error degrees of freedom for the DMM tests. Boik (1991) 
proposed an alternative adjustment factor that improves the approximation when Cl has the general Kronecker 
structure given in (2.22). 

Using Wilks’s A criterion, Boik (1988) shows how one may make the ^-adjustment to the associated F- 
statistic. The result applies equally to the other criteria using Rao’s F-approximation. In Program l.sas, we use , 



the Bartlett-Lawley-Hotelling trace criterion. Then 



F = 2(sN+l)T^ /s^(2M+s+ 1) 
- f[j(2A/ + j+1),2(jA + 1)] 



(4.18) 



where 

s= min (vl,p),vl = h4v^,vl=h4Vf, 

M=(\vl~ pl-1) l2,N=iv*-p-\)l2. 

From the program output, £ = 0.7305055 and the £- adjusted MANOVA Table 3 for the MMM tests results. 
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Table 3 



, e - adjusted Analysis, 


Zullo Data 






Hypothesis 


F 


df 


p-value 


Cond/Parall 


65.090 


(4.38, 30.31) 


< 0.0001 


Group X Cond 


0.903 


(4.38, 30.31) 


0.4830 



The p-value for the unadjusted MMM was 0.4862 for the test of parallelism. The comparison using Wilks’s A 
criterion is 0.463 versus 0.4687, for the adjusted MMM and the unadjusted test, respectively. The nominal p- values 
for the two procedures are approximately equal, 0.483 vs. 0.463, for and A, respectively, for the test of 
parallelism. 

It appears that we have three competing strategies for the analysis of a DMM. The DMM, the MMM, and the 
e- adjusted MMM. Given multivariate sphericity, the MMM is of course most powerful. When multivariate 
sphericity does not hold, neither the adjusted MMM or the DMM analysis is most powerful. Boik (1991) 
recommends using the DMM. The choice between the two procedures depends on the ratio of the traces of the 
noncentrality matrices of the associated Wishart distributions which is seldom known in practice. 

Extended Linear Hypotheses - Example 

Krishnaiah et al. (1980) developed a “Roots” program to test extended linear hypotheses using three test criteria: 
7^ , Roy’s largest root and T^. They consider hypotheses using only a single matrix G. We illustrate their 
procedure using Roy’s largest root test and using PROC IML and data from Timm (1975, p. 454) shown in Table 
4 involving three groups and repeated measurements data. 
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Table 4 

Repeated Measurements. Timm (1975. d. 4541 





Subjects 


1 


Conditions 
2 3 




1 


2 


4 


7 




2 


2 


6 


10 


Drug group 1 


3 


3 


7 


10 




4 


7 


9 


11 




5 


6 


9 


12 


Means 




4 


7 


10 




1 


5 


6 


10 




2 


4 


5 


10 


Drug group 2 


3 


7 


8 


11 




4 


8 


9 


11 




5 


11 


12 


13 


Means 




7 


8 


11 




1 


3 


4 


7 




2 


3 


6 


9 


Drug group 3 


3 


4 


7 


9 




4 


8 


8 


10 




5 


7 


10 


10 


Means 




5 


7 


9 


Grand mean 




5.333 


7.333 


10 



Using PROC GLM we test the standard MANOVA test for groups, conditions, and parallelism (group x 
condition interaction) using Program 2.sas in Appendix B. Finding the differences in groups to be significant, 
suppose we are interested in the extended linear hypothesis y = 0 where y is given (3.9). Using the PROC IML, 
we verify the calculations of the overall test and estimate 

y = rr(G0) = rr(GC„BA„) = Tr(A„GC„B) 

= Tr(G*B) = -2.5 

with the matrices G, G* , and B defined: 



G = 


^ 1 .5> 
-.5 .5 


,G* = 


1 -.5 -.5' 
-.5 1 -.5 


and B = 


^4 7 10"^ 
7 8 11 




,-.5 1, 




-.5 -.5 1 




,5 7 9, 
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Solving IH - 5 E;' 1=0 with H = GW„G', W„ = C„(X'X)''C; and E„=A;EA„ = 

A' Y'(I-X(X'X)“'X')YA„, we solve the characteristic equation using the IML routine EIGVAL for a symmetric 
matrix. Because of rounding, the last characteristic root is essentially zero and, 5] = 6.881 and ^2 = 2.119. By 
(3.4), 

^Lcc = 7’r(GW„GEj=z5,. 

I 

for the Bartlett-Lawley-Hotelling trace criterion. The corresponding value of a for largest root criterion is 

I 

as shown by Krishnaiah et al. (1980). 

Since ij/ = -25, the extended trace and root statistics for testing y/ = 0 are: 

I^I/a^race = 0.8333 
\y/\/a^, = 1.3320. 

Evaluating (2.13) with a = 0.05 for the trace criterion and a similar result for Roy’s largest root test, the critical 
values for the criteria are 1.331991 and 0.9891365, respectively. As shown in the program output, approximate 
confidence intervals for y/ are (-53.77, 48.78) and (-6.50, -1.50) indicating nonsignificance of the contrast y/ for 
both criteria. 

We next illustrate testing (3.10) using the PROC IML by evaluating the supremum in (3.5) as outlined in (3.6). 
With Ti = 7>(G;0) and ty = Tr(Gi Gy E^), we have in the Output that 

t’=(-3,-1, 1, 2) 



^ 29.6 26 -13 -7 ^ 

26 27.2 -13.6 -7.6 

-13 -13.6 27.2 15.2 

, 7 -7.6 15.2 10.4; 

so that t' T“’t = 2.2081 151. Comparing this with the U(G) critical value (1.3319921)^ = 1.774203, we see that the 
test is significant. To find confidence intervals for contrasts involving y/i = 7>(G,.0), one would find a 1-a 
simultaneous confidence interval for each contrast as illustrated above. Finally, 95% confidence intervals are 
constructed for \|f in (3.11) following the test of parallelism, and for y/ in (3.12) following the test of equal 
condition vectors. Program 2.sas. 

5. An extended linear MMM hypothesis test 

As the within subject design in a repeated measurements experiment becomes more complex, we often find that 
the multivariate sphericity condition is not satisfied. For complex designs, the more general form of the hypothesis 
in (2.3) becomes 

H:0 = CB(L(g)A) = O (5.1) 
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where is an orthogonal contrast matrix of rank p,L'L = Given (5.1), the covariance structure for (2.1) is: 

f2 = (L® A)'I(L®A) (5.2) 

and 0 = CB(L® A). Using the DMM, one may test (5.1) or (3.2) for arbitrary 

We can also test (5.1) or (3.2) given multivariate sphericity. Suppose, however, that does not satisfy 
the multivariate sphericity condition, but has general Kronecker structure: 

= (5.3) 



where L^(px p) and = A'ZA(wxw) are arbitrary positive definite matrices, a structure commonly found in 
three-mode factor analysis models, Rentier and Lee (1978). 

To test the null hypothesis that has the structure given in (5.3) versus the alternative that the rows 
of Y have a general structure, under multivariate normality, we obtain a likelihood ratio test using the normal 



likelihood, Timm (1975, p. 558). The likelihood ratio statistic for testing the covariance structure is given by 



(5.4) 



whereZ^ and Z„ are the maximum likelihood (ML) estimates of and Z„ under the null hypothesis, 

N ^ N 

S = 2 (y,. - y)(y|- - y)' / ^ is the ML estimate of the covariance matrice and y = Z yj / ^. 

i=l 1=1 

To test for the structure, the statistic ~2\nX is asymptotically distributed as a chi-square distribution with 
V degrees of freedom where v = p^(p^ + l)/2-[p(p + l) + ^(^ + l) + l]/2 = (p-l)(^-l)[(p + l)(^ + 1)+ 1]/ 2. 

However, to solve the likelihood equations to obtain Z^ andZ„ involves an iterative process as outlined by 
Krishnaiah and Lee (1980), Boik (1991),and Naik and Rao (1996). Naik and Rao provide a computer program 
using the SAS IML procedure to obtain the ML estimates. 

Given that satisfies (5.3) and Z^ and Z„ are estimated by Z^ and Z„ , we may develop a test of the 
extended linear hypothesis 

H:v^ = rr(G0) = O. (5.5) 



Following Mudholkar et al. (1974), a test statistic for testing (5.5) is 

= {(/(G)}^ = [rKG©)]^ / rr[GC(X'X)"' C'G(i, ® i„ )] (5.6) 



since 

w ~ N(y/,al = rr[GC(X'X)’'C'G(i, ®ij] (5.7) 

for fixed G. The statistic converges to a chi-square distribution with one degree of freedom. 

The statistic in (5.7) is an alternative to Boik’s e- adjusted multivariate test procedure for the more general 
hypothesis given in (5.1) when multivariate sphericity is not satisfied. While we have provided a asymptotic test of 
(5.1) given (5.3), the more difficult problem is the estimation of Z^ and Z„, a solution to the likelihood in (5.4) 
over all positive definite matrices Z^ and Z^. Naik and Rao (1996) have developed an alternative Satterthwaite 
type approximate for the MANOVA model when multivariate sphericity is not satisfied. 
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Appendix A 



/* Program l.sas */ 

/* Example from Timm(1980) and Boik(1988,1991) - Zullo Dental data */ 
options LS=78 ps=60 nodate nonumber; 

title 'Output: Double Multivariate Linear Model and Multivariate Mixed Model’; 
data dmlm; 
infile 'mmm.dat'; 
input group yl - y9; 
proc print data=dmlm; 
proc glm; 
class group; 

model yl - y9 = group/nouni; 
means group; 

/* Multivariate test of group differences for mean vectors */ 
manova h=group /printh printe; 
proc glm; 
class group; 

model y 1 - y9 = group / nouni; 

/* Multivariate test of Parallelism */ 
contrast Tarallel' 
group 1-1; 

manova m=(.7071 0 -.7071 0 0 0 0 0 0, 

0 0 0. -.408 .816 -.408 0 0 0, 

00 0 0 00.70710 -.7071, 

-.408 .816-.408 000 0 0 0, 

0 0 0.70710 -.7071 0 0 0, 

00 0 000 -.408 .816 -.408) 

prefix = pari/ printh printe; 

proc glm; 
class group; 

model y 1 - y9 = group /noint nouni; 

/* Multivariate test of Conditions as vectors */ 
contrast 'Mult Cond' group 1 0, 
group 0 1; 

manova m=(l -10 0 0 0 0 0 0, 

0 1 -1 000000 , 

000 1-10000, 

0000 1 - 1000 , 

000000 1 -1 0 , 

00000001-1) prefix = diff^printh printe; 

proc glm; 
class group; 

model y 1- y9 = group/noint nouni; 

/* Multivariate test of Conditions given Parallelism */ 
contrast ’CondlParl’ group .5 .5; 
manova m=(.707 1 0 -.707 1 00 0 0 0 0, 

000 -.408 .816 0 -.408 0 0 0, 

000000.7071 0 -.7071, 

-.408 .816 -.408 000000, 

000-70710 -.7071 00 0, 

0 0 0 0 0 0 -.408 .816 -.408) prefix=cond/printh printe; 

proc glm; 
class group; 

model y 1 - y9 = group/noint nouni; 
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contrast 'GrouplParr group 1-1; 



manova m=(.577 .577 .577 0 00 0 0 0, 

0 0 0 .577 .577 .577 0 0 0, 

0 0 0 000 .577 .577 .577) prefix=Ovall/printh printe; 

/* Multivariate Mixed Model Analysis */ 
data mix; 

infile 'mixed.dat'; 
input group subj cond yl y2 y3; 
proc print data=mix; 
proc glm; 

class group subj cond; 

model yl - y3 = group subj(group) cond cond*group; 
random subj (group); 

contrast ’Group’ group 1 -l/e=subj (group); 
manova h = cond group*cond/printh printe; 

/* Test for Multivariate Spericity and calculation of Epsilon for MMM */ 
proc iml; 

print ’Test of Multivariate Sphericity Using Chi-Square and Adjusted Chi- 
Square Statistics'; 

e={ 9.6944 7.3056 -6.7972 -4.4264 -0.6736 3.7255, 

7.3056 8.8889 -4.4583-3.1915 -3.2396 2.9268, 

-6.7972-4.4583 18.6156 2.5772 0.8837-10.1363, 

-4.4264-3.1915 2.5772 5.3981 1.4259 -1.8546, 

-0.6736-3.2396 0.8837 1.4259 18.3704 -.7769, 

3.7255 2.9268-10.1363-1.8546-0.7769 6.1274}; 
print e; 
n=18; 
p=3; 
t=3; 
k=2; 
u=6; 
q=u/p; 
nu_e=n-k; 
nu_h=l; 

el l=e[l:3,l:3];print ell; 

e22=e[4:6,4:6];print e22; 

dn=(ell+e22)/2; 

b=eigval(dn); print b; 

a=eigval(e); print a; 

b=log(b); 

a=log(a); 

chL2=n#(q#sum(b)-sum(a)); 
df=p#(q- l)#(p#q+p+l)/2; 
pvalue=l-probchi(chi_2,df); 
print chi_2 df pvalue; 
cl=p/(12#q#nu^e#df); 

rho= 1 -c 1 #(2#p##2#(q##4- 1 )+3#p#(q##3- 1 )-(q##2- 1 )) ; 
ro_chi_2=(rho#nu_e/n)#chi_2; print rho; 
c2=l/(2#rho##2); 

c3=((p#q- 1 )#p#q#(p#q+ 1 )#(p#q+2))/(24#nu_e##2); 
c4=((p- l)#p#(p+ 1 )#(p+2))/(24#q##2#nu_e##2); 
c5=df#(l-rho)##2/2; 
omega=c2#(c3-c4-c5); print omega; 
p 1 = 1 -pr obchi (ro_c hi_2 ,df) ; 




20 



22 



p2=l-probchi(ro_chL2,df+4); 
cp value=( 1 -omega)#p 1 +omega#p2; 
print ro_chi_2 cpvalue; 
s=e/nu_e; 

sll=s[l:3,l:3];sl2=s[l:3,4:6]; 
s21=s[4:6,l:3]; s22=s[4:6,4:6]; 

enum=trace((sll+s22)*(sll+s22))+(trace(sll+s22))##2; 
eden=q#( trace(sl l)##2+trace(sl l*sl l)+trace(sl2)##2+trace(sl2*sl2)+ 
trace(s21)##2+trace(s21*s21)+trace(s22)##2+trace(s22*s22)); 
epsilon=enum/eden; 

nu_h=nu_h#q; nu_e=nu_e#q; sO=min(nu_h,p); 

Mnu_h=nu_h#epsilon; Mnu_e=nu_e#epsilon; msO=min(mnu_h,p); 
mO=(abs(mnu_h-p)- 1)/2; nO=(mnu_e-p- 1)/2; 
denom=sO##2#(2#mO+msO+ 1 ); numer=2#(msO#nO+ 1 ); 
dfl=msO#(2#mO+msO+l);df2=2#(msO#nO+l); 

print ’Epsison adjusted F-Statistics for cond and group X cond MMM tests'; 
print epsilon; 

f_cond = df2#13.75139851/(ms0#dfl); 
f_gXc =df2#0.19070696/(ms0#dfl); 
print f_cond f_gXc df 1 df2; 
p_cond= 1 -probf(f_cond»df 1 ,df2); 
p_gXc=l-probf(f_gXc,dfl,df2); 

print 'Epsilon adjusted pvalues for MMM tests using T0**2 Criterion’; 
print p_cond p_gXc; 
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Appendix B 



* Program 2.sas */ 

/* Data from Timm(1975, page 454) */ 
options ls=78 ps = 60 nodate nonumber; 
title 'Output; Extended Linear Hypotheses' ; 
data timm; 
infile 'exlin.dat'; 
input group y 1 y2 y3 xl x2 x3; 
proc print data=timm; 
proc glm; 
class group; 

model yl-y3 = group/nouni; 
means group; 

/* Multivariate test of Groups */ 
manova h=group/printh printe; 
proc iml; 
use Timm; 
a={xl x2 x3}; 
b=iyl y2 y3}; 
read all var a into x; 
read all var b into y; 
beta=inv(x' *x)*x' *y; 
print beta; 
n=nrow(y); 
p=ncol(y); 
k=ncol(x); 

nu_h=2; u=3; nu_e=n-k; sO=min(nu_h»u); r=max(nu_h»u); alpha=.05; 

denr=(nu_e-r+nu_h); 

roy_2=(r/denr)*finv(l-alpha,r,denr); 

rvalue=sqrt(roy_2); 

mO=(abs(nu_h-u)- 1)/2; n0=(nu_e-u- 1)/2; 

num=s0**2*(2*m0+s0+l); dent=2*(s0*n0+l);df=s0*(2*m0+s0+l); 

t0_2=(num/dent)*finv(l-alpha,df,dent); 

tovalue=sqrt(t0_2); 

print sO mO nO; 

e=(y' *y-y' *x*beta); 

co={l -10,0 1-1}; 

ao=i(3); 

eo=ao'*e*ao; 

bo=co*beta*ao; 

wo=co* in v(x' *x) *co' ; 

ho=bo'*inv(wo)*bo; 

print,"Overall Error Matrix", eo ,"Overall Hypothesis Test Matrix",ho; 

/* c'c=eo where c is upper triangle Cholesky matrix */ 

c=root(eo); 

f=inv(c')*ho*inv(c); 

eig=Eig val(round(f ,.000 1 )) ; 

vec=i n v(c) *eig vec(round(f, .000 1 )) ; 

print,"Eigenvalues & Eigenvectors of Overall Test of Ho (Groups)", eig vec; 
/* Extended Linear Hypothesis following Overall Group Test */ 
m={1.5,-.5 .5, -.5-1}; 
g=ao*m*co; 

print, "Extended Linear Hypothesis Test Matrix",m g; 
psi=m*bo; 





psi_hat=trace(psi); 

tr_psi=abs(psi_hat); 

h=m*wo*m'; 

eo=inv (eo); 

c=root(eo); 

f=inv(c')*h*inv(c); 

xeig=Eig val (round(f, .000 1 )) ; 

print, "Eigenvalues of Extended Linear Hypothesis”, xeig; 
to_2=tr_psi/sqrt(sum(xeig)); print, "Extended To**2 Statistic”, to_2; 
print, "Extended To**2 Critical Value", tovalue; 

root=tr_psi/sum(ssq(xeig)); print, "Extended Largest Root Statistic", root; 

print, "Extended Largest Root Critical Value”, rvalue; 

print psi_hat alpha; 

ru=psi_hat+rvalue*sum(ssq(xeig)); 

rl=psi_hat-rvalue*sum(ssq(xeig)); 

vu=psLhat+tovalue*sqrt(sum(xeig)); 

vl=psi_hat- tovalue* sqrt(sum(xeig)); 

print ’Approximate Simultaneous Confidence Intervals'; 

print 'Contrast Significant if interval does not contain zero'; 

print 'Extended Root interval: (’rl ’,' ru ')'; 

print 'Extended Trace interval: ('vl ',' vu ')'; 

/* Multiple Extended Linear Hypothesis using To* *2 */ 

ml={ 1 0,0 0,0 0}; m2={0 0,1 0,0 0}; m3={0 0,0 1,0 0}; m4={0 0,0 0,0 1 }; 

print,”Multiple Extended Linear Hypothesis Test Matrices", ml,m2,m3,m4; 

g l=ao*m 1 *co; g2=ao*m2*co; g3=ao*m3*co; g4=ao*m4*co; 

tl=trace(ml*bo); t2=trace(m2*bo); t3=trace(m3*bo); t4=trace(m4*bo); 

tau=tl//t2//t3//t4; 

til =trace(m 1 * wo*m 1 ' *eo); 

t2 1 =trace(m2* wo*m 1 ' *eo); t22=trace(m2* wo*m2' *eo); 

t31=trace(m3*wo*mr*eo); t32=trace(m3*wo*m2'*eo); t33=trace(m3*wo*m3'*eo); 

t41=trace(m4*wo*mr*eo); t42=trace(m4*wo*m2'*eo); t43=trace(m4*wo*m3'*eo); 

t44=trace(m4*wo*m4'*eo); 

rl=tllllt21llt31llt41; 

r2=t21llt22llt32llt42; 

r3=t31llt321lt33llt43; 

r4=t41llt42llt43llt44; 

t=rl//r2//r3//r4; 

print tau,t; 

to_4=tau'*inv(t)*tau; 

print, "Extended Linear Hypothesis Criterion To**2 Squared", to_4; 
print, "Extended To**2 Critical Value", t0_2; 

/* Multivariate test of Parallelism */ 
data timm; 
infile 'exlin.dat'; 
input group yl y2 y3 xl x2 x3; 
proc glm; 
class group; 

model yl-y3 = group/nouni; 
manova h = group m = ( 1 -1 0, 

0 1-1) prefix = dif^printe printh; 

proc iml; 
use Timm; 
a={xl x2 x3}; 
b={yl y2 y3}; 
read all var a into x; 
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read all var b into y; 

beta=inv(x'*x)*x'*y; 

n=nrow(y); 

p=ncol(y); 

k=ncol(x); 

nu_h=2; u=2; nu_e=n-k; sO=min(nu_h,u); r=max(nu_h,u); alpha=.05; 

denr=(nu_e-r+nu_h); 

roy_2=(r/denr)*finv( 1 -alpha,r,denr); 

rvalue=sqrt(roy_2) ; 

m0=(abs(nu_h-u)-l)/2; n0=(nu_e-u-l)/2; 

num=s0**2*(2*m0+s0+l); dent=2*(s0*n0+l); df=s0*(2*m0+s0+l); 

t0_2=(num/dent)*finv(l -alpha, df, dent); 

tovalue=sqrt(t0_2); 

print sO mO nO; 

e=(y'*y-y'*x*beta); 

co={l -1 0,0 1 -1}; 

ao={ 10,-1 1,0-1}; 

eo=ao'*e*ao; 

bo=co*beta *ao; 

wo=co*inv(x' *x)*co^ ; 

ho=bo'*inv(wo)*bo; 

c=root(eo); 

f=inv(c')*ho*inv(c); 

eig=eigval(round(f,.0001)); 

vec=inv(c)*eigvec(round(f, .000 1 )); 

print, "Eigenvalues & Eigenvectors of Overall test of Ho (Parallelism)", eig vec; 
/* Extended Linear Hypothesis following overall Parallelism test */ 
m={01,l 0); 
g=ao*m*co; 

print, "Extended Linear Hypothesis Test Matrix", m g; 

psi=m*bo; 

psi_hat=trace(psi); 

tr_psi=abs(psi_hat); 

h=m*wo*m'; 

eo=inv(eo); 

c=root(eo); 

f=inv(c^)*h*inv(c); 

xeig=eig val(round(f,.000 1 )) ; 

print, "Eigenvalues of Extended Linear Hypothesis", xeig; 
to_2=tr_psi/sqrt(sum(xeig)); print, "Extended To**2 Statistic", to_2; 
print,"Extended To**2 Critical Value", tovalue; 

root=tr_psi/sum(ssq(xeig)); print, "Extended Largest Root Statistic", root; 

print, "Extended Largest Root Critical Value", rvalue; 

print psi_hat alpha; 

ru=psLhat+rvalue*sum(ssq(xeig)); 

rl=psi_hat-rvalue*sum(ssq(xeig)); 

vu=psi_hat+tovalue*sqrt(sum(xeig)); 

vl=psLhat-tovalue*sqrt(sum(xeig)); 

print 'Approximate Simultaneous Confidence Intervals'; 

print 'Contrast Significant if interval does not contain zero'; 

print 'Extended Root Interval: ('rl ru ')'; 

print 'Extended Trace Interval: ('vl vu ')'; 

/* Multivariate test of Conditions as vectors */ 
data timm; 
infile 'exlin.dat'; 
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input group yl y2 y3 xl x2 x3; 
proc glm; 
class group; 

model yl-y3 = group/noint nouni; 
contrast 'Mult Cond' group 100, 
group 0 10, 
group 00 1; 

manova m=(l -10, 

0 1-1) prefix = diff/ printe printh; 

proc iml; 
use timm; 
a={xl x2 x3}; 
b={yl y2 y3}; 
read all var a into x; 
read all var b into y; 
beta=inv(x' *x)*x' *y ; 
n=nrow(y); 
p=ncol(y); 
k=ncol(x); 

nu_h=3; u=2; nu_e=n-k; s0=min(nu_h,u); r=max(nu_h,u); alpha=.05; 

denr=(nu_e-r+nu_h); 

roy_2=(r/denr) *fin v( 1 -alpha,r,denr); 

rvalue=sqrt(roy_2); 

m0=(abs(nu_h-u)- 1 )/2; n0=(nu_e-u- 1)/2; 

num=s0**2*(2*m0+s0+l); dent=2*(s0*n0+l); df=s0*(2*m0+s0+l); 

t0_2=(num/dent)*finv(l-alpha,df,dent); 

to value=sqrt(t0_2) ; 

print sO mO nO; 

e=(y' *y-y' *x*beta); 

co=i(3); 

ao={ 10,-1 1,0-1}; 

eo=ao'*e*ao; 

bo=co*beta*ao; 

wo=co*inv(x' *x)*co' ; 

ho=bo'*inv(wo)*bo; 

c=root(eo); 

f=inv(c')*ho*inv(c); 

eig=eigval(round(f,.0001 )); 

vec=inv(c)*eigvec(round(f,.0001)); 

print, "Eigenvalues & Eigenvectors of Overall test of Ho (Conditions)", eig vec; 

m={101,01 1); 

g=ao*m*co; 

print, "Extended Linear Hypothesis Test Matrix", m g; 

psi=m*bo; 

psi_hat=trace(psi) ; 

tr_psi=abs(psLhat); 

h=m*wo*m'; 

eo=inv(eo); 

c=root(eo); 

f=inv(c')*h*inv(c); 

xeig=eig val(round(f, .000 1 )) ; 

print, "Eigenvalues of Extended Linear Hypothesis", xeig; 
to_2=tr_psi/sqrt(sum(xeig)); print, "Extended To**2 Statistic", to_2; 
print, "Extended T0**2 Critical Value", tovalue; 

root= tr_psi/sum(ssq(xeig)); print, "Extended Largest Root Statistic", root; 
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print, "Extended Largest Root Critical Value", rvalue; 

print psLhat alpha; 

ru=psLhat+rvalue*sum(ssq(xeig)); 

rl=psLhat-rvalue*sum(ssq(xeig)); 

vu=psLhat+tovalue*sqrt(sum(xeig)); 

vl=psLhat-tovalue*sqrt(sum(xeig)); 

print ’Approximate Simultaneous Confidence Interveds'; 

print ’Contrast Significant if interval does not contain zero'; 

print ’Extended Root Interval: (’rl ru ’)’; 

print ’Extended Trace Interval: (’vl ’,’ vu ’)’; 
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